
;;cnum = '4594'
cnum = '1724'
ellcol = 'Powder Blue'
ellth = 5.

orig = readfits('./local/ubc_glimpse_proc/nonsmooth'+cnum+'.fits',ohdr,/sil)
smoo = readfits('./local/ubc_glimpse_proc/smooth30arc'+cnum+'.fits',shdr,/sil)
imir = readfits('./local/ubc_glimpse_proc/Imir'+cnum+'.fits',ihdr,/sil)

restore,'irdc_full_data.sav',/ver
s = read_bgps_csv(/ver)
i = WHERE(s.cnum EQ long(cnum),ni)

;; Get full GLIMPSE frame
glon = string( round(s[i].glon_peak / 3.) * 3., format="(I03)")
;gimg = readfits('./local/glimpse_v3/GLM_'+glon+'00+0000_mosaic_I4.fits',ghdr)
gimg = readfits('./local/glimpse_v3/starsub.'+glon+'00.I4.fits',ghdr)
plot_wcs_axes,ohdr,10,lonarr=xc,latarr=yc



gcirc,2,s[i].glon_peak,s[i].glat_peak,pf_irdc.l,pf_irdc.b,dis

jj = WHERE(dis LE 105., njj)
print,njj



myps,'./irdc_dist_model/analysis_plots/ubc_example_fig.eps',xsize=11,ysize=10

multiplot_xm, [2,2], mtitle='BGPS #'+cnum,gap=0.03,mtitoffset=0, xmargin=1

pr = set_plot_range(alog(smoo))

plotmap,gimg,ghdr,charsize=0.8,ct=0,range=pr,xrange=[max(xc),min(xc)],$
        yrange=[min(yc),max(yc)],/log

FOR k=0L, njj-1 DO BEGIN
   j = jj[k]
   tvellipse,pf_irdc[j].min/3600.,pf_irdc[j].maj/3600.,$
             pf_irdc[j].l,pf_irdc[j].b,$
             -pf_irdc[j].pa_irdc,/data,thick=ellth,color=cgColor(ellcol)
ENDFOR

wtmask  = aperture_wtmask(ghdr, 40., s[i].glon_peak, s[i].glat_peak)
plot_wcs_axes,ghdr,10,lonarr=xc,latarr=yc
contour,wtmask,xc,yc,levels=[0.5],thick=5,color=cgColor('Deep Pink'),/overplot

al_legend,/top,/right,['(a)'],textcolor='Spring Green',charsize=1.5,box=0

multiplot


plotmap,orig,ohdr,charsize=0.8,ct=0,range=pr,/log,min_dpi=100


FOR k=0L, njj-1 DO BEGIN
   j = jj[k]
   tvellipse,pf_irdc[j].min/3600.,pf_irdc[j].maj/3600.,$
             pf_irdc[j].l,pf_irdc[j].b,$
             -pf_irdc[j].pa_irdc,/data,thick=ellth,color=cgColor(ellcol)
ENDFOR

wtmask  = aperture_wtmask(ohdr, 40., s[i].glon_peak, s[i].glat_peak)
plot_wcs_axes,ohdr,10,lonarr=xc,latarr=yc
contour,wtmask,xc,yc,levels=[0.5],thick=5,color=cgColor('Deep Pink'),/overplot


al_legend,/top,/right,['(b)'],textcolor='Spring Green',charsize=1.5,box=0

multiplot

plotmap,smoo,shdr,charsize=0.8,ct=0,range=pr,/log,min_dpi=100

FOR k=0L, njj-1 DO BEGIN
   j = jj[k]
   tvellipse,pf_irdc[j].min/3600.,pf_irdc[j].maj/3600.,$
             pf_irdc[j].l,pf_irdc[j].b,$
             -pf_irdc[j].pa_irdc,/data,thick=ellth,color=cgColor(ellcol)
ENDFOR

wtmask  = aperture_wtmask(shdr, 40., s[i].glon_peak, s[i].glat_peak)
plot_wcs_axes,shdr,10,lonarr=xc,latarr=yc
contour,wtmask,xc,yc,levels=[0.5],thick=5,color=cgColor('Deep Pink'),/overplot

al_legend,/top,/right,['(c)'],textcolor='Spring Green',charsize=1.5,box=0

multiplot

plotmap,imir,ihdr,charsize=0.8,ct=0,range=pr,/log


cgLoadct,0,/silent
cgcolorbar,range=exp(pr),/vertical,/right,position=[0.92,0.09,0.94,0.94],/ylog,$
           format="(F0.1)",title='GLIMPSE Intensity [MJy/sr]',charsize=0.8

FOR k=0L, njj-1 DO BEGIN
   j = jj[k]
   tvellipse,pf_irdc[j].min/3600.,pf_irdc[j].maj/3600.,$
             pf_irdc[j].l,pf_irdc[j].b,$
             -pf_irdc[j].pa_irdc,/data,thick=ellth,color=cgColor(ellcol)
ENDFOR

wtmask  = aperture_wtmask(ihdr, 40., s[i].glon_peak, s[i].glat_peak)
plot_wcs_axes,ihdr,10,lonarr=xc,latarr=yc
contour,wtmask,xc,yc,levels=[0.5],thick=5,color=cgColor('Deep Pink'),/overplot

al_legend,/top,/right,['(d)'],textcolor='Spring Green',charsize=1.5,box=0

multiplot,/reset
multiplot,/default
multiplot,/reset

myps,/done



END
